### ### ### ### ### ### ###  Replication script ERG models ### #### ### ### ### ### ### 

### load data
load("XXXX/Replication_Data_Models.RData")

### load library
library("statnet")
library("stargazer")

########################################################################################
####################################### Models #########################################
sink('ERGM-output.txt')
seed <- 12345
set.seed(seed)

# best practice model

final_practice = ergm(nw.practice ~ edges + gwesp(.3, fixed = TRUE) + 
                         nodematch("welfare")+
                         edgecov(geo_mat) +
                         nodecov("staff") +
                         nodecov("receiving") +
                         nodecov("sending") +
                         nodecov("capacity"),
                       control = control.ergm(MCMC.samplesize=50000,MCMC.interval=1000),
                       verbose=T)

summary(final_practice)

# advice model

final_advice = ergm(nw.advice_sym ~ edges + gwesp(.3, fixed = TRUE) + # gwdsp(.3, fixed = TRUE) +
                      nodematch("welfare")+
                      edgecov(geo_mat) +
                      nodecov("staff") +
                      nodecov("receiving") +
                      nodecov("sending") +
                      nodecov("capacity"),
                    control = control.ergm(MCMC.samplesize=50000,MCMC.interval=1000),
                    verbose=T)

summary(final_advice)


# information model

final_info = ergm(nw.information ~ edges + gwesp(.3, fixed = TRUE) + # gwdsp(.3, fixed = TRUE) +
                    nodematch("welfare")+
                    edgecov(geo_mat) +
                    nodecov("staff") +
                    nodecov("receiving") +
                    nodecov("sending") +
                    nodecov("capacity"),
                  control = control.ergm(MCMC.samplesize=50000,MCMC.interval=1000),
                  verbose=T)

summary(final_info)

# problem solving model

final_problem = ergm(nw.problem ~ edges + gwesp(.3, fixed = TRUE) + 
                       nodematch("welfare")+
                       edgecov(geo_mat) +
                       nodecov("staff") +
                       nodecov("receiving") +
                       nodecov("sending") +
                       nodecov("capacity"),
                     control = control.ergm(MCMC.samplesize=50000,MCMC.interval=1000),
                     verbose=T)

summary(final_problem)


### Table 1
stargazer(final_practice, final_advice, final_info, final_problem, type="text")


### Goodness of fit diagnostics
gof_practice <- gof(final_practice ~ degree + esp + distance)
par(mfrow = c(2,2))
plot(gof_practice, main= "Goodness-of-fit diagnostics: Best practices")
gof_advice <- gof(final_advice ~ degree + esp + distance)
par(mfrow = c(2,2))
plot(gof_advice, main= "Goodness-of-fit diagnostics: Advice")
gof_info <- gof(final_info ~ degree + esp + distance)
par(mfrow = c(2,2))
plot(gof_info, main= "Goodness-of-fit diagnostics: Information")
gof_problem <- gof(final_problem ~ degree + esp + distance)
par(mfrow = c(2,2))
plot(gof_problem, main= "Goodness-of-fit diagnostics: Problem solving")

sink()
